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Abstract 

Phase space method provides a novel way for deducing qualitative features 
of nonlinear differential equations without actually solving them. The method 
{^j^l is applied here for analyzing stability of circular orbits of test particles in var- 

ious physically interesting environments. The approach is shown to work in a 
revealing way in Schwarzschild spacetime. All relevant conclusions about circu- 
^ ■ lar orbits in the Schwarzschild-de Sitter spacetime are shown to be remarkably 

0^ ', encoded in a single parameter. The analysis in the rotating Kerr black hole 

. readily exposes information as to how stability depends on the ratio of source 

rotation to particle angular momentum. As a wider application, it is exempli- 
fied how the analysis reveals useful information when applied to motion in a 

op ' refractive medium, for instance, that of optical black holes. 
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1. Introduction 



It is well known that in Einstein's theory of General Relativity (GR) motion 
' in a gravitational field is described by a system of dynamical equations. The 

nonlinear ordinary differential equation for the path is obtained by eliminating 
the affine parameter from that system of equations. An important class of solu- 
tions of the path equation is formed by circular trajectories. The issue of their 
stability is particularly important in confinement problems and/or in accretion 
phenomenon in astrophysics. However, conventional analysis of stability of or- 
bits essentially deals with dynamical equations involving the affine parameter 
and a potential function. On the other hand, potential functions may not always 
be immediately evident, for instance, in a simulated environment like moving 
refractive dielectrics whereas path equations in them could follow directly from 
Fermat's principle or Hamilton- Jacobi equation. Therefore, a natural query is 
to ask if information about stability of circular orbits can be obtained from the 
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geometrical path equation alone. Our aim here is to rigorously demonstrate 
that it is indeed possible via phase space analysis of autonomous systems cor- 
responding to various physically interesting environments. To our knowledge, 
such a useful application in gravitational physics seems yet unavailable in the 
literature. 

For our analysis, we shall not require any other information beyond the 
path equation. The conserved angular momentum appearing in it will allow 
us to connect the phase space results with actual kinematics. An equilibrium 
state corresponds to a constant solution of a differential equation describing a 
mechanical system and conversely. Constant solution means that velocity x and 
acceleration x be simultaneously zero. The concept of stability of an equilibrium 
state is borrowed from the familiar example of motion of a pendulum about the 
equilibrimxi point x = and a; = where x is angle with the vertical. The 
motion is stable because a small displacement from the equilibrium position 
will lead to only small oscillations of the bob about that position in a vertical 
plane. Different closed paths on the phase space about a stable equilibrium 
point correspond to real oscillations with different periods. If the paths around 
an equilibrium point is such that a small displacement from the equilibrium state 
takes the system far away from that point, it is called an unstable equilibrium 
point. For instance, the equilibrium point x = tt and i = is an unstable saddle 
in the pendulum motion. There is no closed path around such equilibrium 
point although some paths may pass through it depending on the values of 
the parameter. Open paths not passing through an equilibrium point represent 
whirling motion of the pendulum [1]. 

In this paper, we shall portray path equations for different solutions of GR 
as plane autonomous dynamical systems and study them through phase space 
and/or Hamiltonian analysis. Of particular interest is the treatment in cosmo- 
logical scenario and in the Kerr metric. As a further example, the method is 
applied to light propagation in a medium relevant to optical black holes and 
interesting information obtained. 

The article is intended for theoretical physicists in general and practising 
relativists in particular. The contents are organized as follows. In Sec. 2, we 
briefly describe the path equation needed for our purpose. Sec.3 treats the path 
equation as an autonomous system without the cosmological constant A = 0. 
Sec. 4 treats the autonomous system as a Hamiltonian system. Sec.5 deals with 
the case A ^ 0. In Sec. 6, we discuss circular motion in Kerr spacetime. Sec.7 
shows an example applying the method to light motion in a refractive medium 
relevant to optical black holes. Sec. 8 summarizes the obtained results. We take 
G = Co = 1, unless specifically restored. 

2. Path equation 

A spherically symmetric static solution of the Einstein field equations for 
a mass M is given by the Schwarzschild-de Sitter solution (SdS) in standard 
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coordinates {x°')={t,r,6,(j)): 



dr^ = B{r)dt^ -B-^{r)dr'^ -r'^d.e'^ -r'^sm^edcj)^, (1) 
r 3 

where A > is the cosmological constant. The cosmological constant A ~ 
10~^^cm~^ > is responsible for dark energy that explains the currently ob- 
served accelerated cosmic expansion [2]. The case A < (anti-de Sitter) is 
observationally ruled out and we shall not deal with this case here. The case 
M 7^ 0, A = corresponds to the pure Schwarzschild solution and tr = 2M is 
the horizon radius so that the metric is valid for r > tr- When M = 0, A 
we have pure de Sitter solution which can be reexpressed in such a manner that 
it represents an expanding space devoid of matter. 

When M 7^ 0, A / 0, there occur two horizons provided that < A < 
^crit = g]g2 • The black hole horizon appears at = ^^j^2 ^'^^ ^^T^ ^^'^ 
cosmological horizon at rc = "T^f cos where ^ = cos~^(-\/3AM^). The 
spacetime is dynamic for r < rh and for r > rc- At the critical value A = \crit, 
the two horizons coincide at r-ph = 3M. The static radius Vst is defined as a 
hypersurface where the attraction due to M balances the cosmic repulsion due 
to A and is given by [3] 

r.= (-) . (2) 

All circular orbits are bounded from below by r = rph and from above by the 
static radius r = rst- 

Following conventional stability analyses of circular orbits involving a po- 
tential function V{r), we obtain the expression 

(fV _ /2M\ /^XR^-5XMR^ + R-QM 

\r=R- - y-^j y i?-3M 

In the case A = 0, it follows from Eq.(3) that stable circular orbits may exist 
only at radu R > 6M. At i? = 6Af , \r=R= which indicates that it is a 
point of inflection. When R +00 and 3M+, we have 0— and +00 

respectively. The first limit indicates stability in the asymptotic region while 
the divergent second limit indicates instability. 

Defining m = i, we have the path equation on the equatorial plane 9 = tt/2 
for a particle as 

X = a + hx -\- cx^ -\- dx~^ (4) 
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where an overdot denotes differentiation with respect to (j) and 



X = 


u = - 




r 




M 


a = 






b = 


-1 


c = 


3M 


d = 


A 






h = 


2# 





(5) 

(6) 

(7) 
(8) 

(9) 

■St. (10) 

The quantity h is the conserved angular momentum. Eq.(4) with the defined 
coefficients is all we need. 

We shall study it in the phase plane {x = u,y = ^) posing it as an au- 
tonomous system as follows 

X = y (11) 
y = a + bx + cx'^ + dx~^ . (12) 

Equilibrium points of the system arc given by 

x = 0,y = 0. (13) 

The first equation gives circular orbits r = R = const, and the second equation 
gives the angular momentum h along that orbit 

,2 MR{1-^R^) ^ ^ 

= ■ (14) 

It follows that = at R = rgt and = oo at R = 3M. Let us now consider 
the case A = corresponding to the Schwarzschild spacetime. 

3. Case I: A = 

The autonomous system 

X = y (15) 
y = a + bx + cx^ (16) 

gives two equilibrium points (a;, t/) at 

In order to have these points located on the real phase plane, we have to assume 
— iac > or rephrasing, = 1 — > 0. We have introduced the 
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shorthand a for notational convenience. Each choice of a gives a corresponding 
equiUbrium point or a value of the radius r. Let us first consider the degenerate 
case a = and study the stabiUty of corresponding radius. 

Case (a): k"^ = b"^ - 4ac = ^ a = |^ 7^ 0. 

The equihbrium point on the phase plane occurs only at (—^,0). Eliminat- 
ing the parameter a, the autonomous system (15), (16) can be reduced to the 
following set of equations 

i = y (17) 
y = ^^{b + 2cxf. (18) 

The differential phase path is given by 

dy {b + 2cxf 



which integrates to 



dx Acy 
1 



(19) 



y^ = ^{b + 2cxf+A (20) 

where A is an arbitrary constant of integration. 

We see that the parameter 6 has the effect of only translation in the variable 
X while c introduces magnification in both x and y. Thus the GR correction 
term c can be regarded as the dominating parameter among a, b and c. By the 
translation 

y' = y (21) 

x' = x+^^ (22) 

the autonomous system (15), (16) further reduces to 

X = y' (23) 
y = cx''^ (24) 

which gives a one-parameter family of phase paths 

3j/'2 = 2cx'^ + C (25) 

on the phase plane {x' ,y') whore C is an arbitrary parameter. The phase paths 
for M = 1 (or c = 3) and different values of C are given in Fig.l. 

The equilibrium point has now been shifted to the origin (0, 0) which gives 

x' = 0^r = 6M^hu = = const. (26) 

dr 

y' = 0^-/.^^?=0. (27) 

^ deb dr ^ ' 
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From the above, we immediately learn the following: The equilibrium point 
corresponds, in the physical (r, </>) plane, to a circular orbit of radius 6M with 
the test particle having a constant cross radial velocity hu and a zero radial 
velocity —h^. From the overall pattern of the phase paths given in Fig.l, we 
see that a small displacement from the equilibrium state can take the system on 
a phase path which leads it far away from the equilibrium state. The dynamical 
condition for this to happen is given by 

= /i^ = 12M2. (28) 

Although the phase path Eq.(25) is independent of a, it applies only to massive 

test particles because the value of h"^ becomes infinity for light (dr = 0). In 
this case, the condition (28) becomes obviously inapplicable. We shall treat this 
case separately in Sec. 4. 

Let us analyze in a little more detail the paths in different quadrants in 
Fig.l. A typical initial state (x'^jUq) on the phase plane is as follows 

x', = xo + ^^=xo--^=6 (29) 




y'o± = yo± = ±f-^^ = ±\r-^ (30) 

where C > —6S^. These equations will allow us to closely examine phase paths 
in the neighborhood of the equilibrium point. We see from Eq.(29) that the 
GR allowed open interval Xq^ = Tq S (2M, +oo) is mapped onto a finite open 
interval x'q G (— g^, 3^) around the equilibrium point (0,0). This interval can 
be subdivided into two parts for 6 or x'q-. 

One part is (5 e [0, 3^) corresponding to 6M > ro > 2M. This interval refers 
to points (xg, ^Q^) on the paths in first quadrant and to points (xq, 2/q_) on the 
paths in the fourth quadrant. To proceed further, let us translate Eqs.(29),(30) 
to the physical (r, </>) plane choosing C = 0: 



6M 

/ ~ rlr 

-h 




d4>) 




0± 



(Note that the cross radial component of velocity can be expressed as = 
hx' and the radial component as ^ = —hy'.) We find the following distinct 
possibilities for paths (25) passing through the equilibrium point (0,0): (i) As 
5 increases from to 3^, we see from Eqs.(29), (30) that both x'^ and y^^ 
increase from the equilibrium point (0, 0) which implies that the phase point in 
the first quadrant moves outward (to the right), as represented by OB in Fig.l. 
Correspondingly, from Eq.(32), we see that u increases with (f) which indicates 
that the radius ro midergocs a decrease in time r from 6M to 2M (as reflected 
in |+= —h^25^ < 0). (ii) As 5 decreases from to 0, we see that both x'q 



6 



and t/Q_ decrease to the equilibrium point (0, 0) which impHes that the phase 
point in the fourth quadrant moves inward (to the left), as represented by AO 
in Fig.l. Correspondingly, u decreases with (f) and the radius tq undergoes an 
increase in time t from 2M to 6M (as reflected in 4^ |_~ +hV2S^ > 0). For 



C ^ 0, paths will not pass through (0, 0) but parts of CDE lying in the first 
and fourth quadrant can be interpreted similarly. 

The other part is (5 G (^el?''^]' ^'^^ which oo > ro > 6M. This represents 
points only on the second and iShird quadrant (where Xq = S <0). In this case, 
to avoid imaginary quantity in Eq.(32), we must choose C ^ 0, that is, we have 
to deal with the full set of Eqs. (29), (30). For different values of C 7^ in Fig.l, 
we see that the phase paths like CDE are not closed around the equilibrium 
point. (These paths are analogous to whirling motion of pendulum). It is clear 
that most of the paths in the interval — < Xq < , when slightly displaced 
from the equilibrium point (0, 0), neither converge to it nor form a center about 
it. In fact, the phase paths resemble exactly those around a cusp [1]. This 
leads us to conclude that the equilibrium point (0, 0) corresponding to radius 
r = 6M is neither stable or nor unstable because of the dynamically degenerate 
condition a = 0. This conclusion will be further supported in Sec. 4. 

Case (b): = - 4ac > 

There are now two distinct equilibrium points occurring at ( "^^^ ,0) and 
( ~2~^ ; 0)- They combine into a single representative point ( ~^^" ,0) where 
a = +k or —k. Under the translation 

y' = y (33) 

(34) 



the autonomous system (15), (16) reduces, after a little algebra, to 

x = y' (35) 

y = ax' + cx''^. (36) 

The equilibrium points in the new {x' , y') phase plane are Pi :(0, 0) and Qi :(— ^, 0). 
The linearized system of equations near Pi :(0,0) is 

X = y' (37) 

y = ax'. (38) 

Comparing it with the general linear system given by 

X = aix' + biy' (39) 

y = cix' + diy' (40) 

we find 

oi = 0,bi = l,ci = a,di = 0,p = ai + di = 0,q = aidi — bici = —a (41) 
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so that the discriminant is 



A = p^ -Aq = 4a. (42) 

Hence the equihbrium point (0, 0) will be a center (stable equilibrium) if a < 
and a saddle point (unstable equilibrium) if a > 0. Such an abrupt change in 
the behavior of the system occurs through a = 0. Therefore, a = can be 
called a bifurcation point. 

The above conclusions are supported by the phase paths following from 
Eqs.(35) and (36), namely, 

= ax'^ + + D (43) 

where D is an arbitrary constant of integration. In the close vicinity of (0, 0) 
such that a;' ~ and x'^ can be neglected, wo get 

y'^ - ax'^ = D (44) 

which represents a family of concentric ellipses for a = —k < (Fig. 2, center, 
stable) and a family of hyperbolas with asymptotes y' = -^^fax' for a = +A; > 
(Fig. 3, saddle point, unstable). Note that a 1, since, for physically realistic 
particle orbits, <C 1. Thus, as a sample, wc have taken a = ±fc = ±0.9, and 
c = 3 (which means we are taking units in which M = 1) and different values 
of the parameter D in Figs. 2, 3. 

Let us sec what conclusions wc can draw in the physical (r, (f)) plane. The 
point (x', 2/')=(0, 0) shows that the equilibrium radii r depend on the value of 
h, and hence of a. These radii follow from Eq.(34) 

x = x — =0=>r=- . (45) 

2c 1 + a ^ ^ 

This means oo > r > 6M if -1 < a < and 6M > r > 2M if < a < 2. 
Thus, from what we have learnt from Eq.(44), we find that circular orbits with 
r > 6M are stable {a has negative sign) while those with r < 6M are unstable 
(a has positive sign). Eq.(33) gives 

/ ^ du ^ dr , du ^ 

while = hu > 0. That is, the radius r is independent of time while the 
cross radial velocity is nonzero. These are exactly what are to be expected 
of circular motions. 

The other equilibrium point {x' ,y')={—^,0) corresponds to 

a a — b 6M 
=x — => r = . (47) 

c 2c 1-a ^ ' 

In the close vicinity of this equilibrium point, we may define 

x" = x' + -. (48) 
c 
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When x" 0, neglecting a;"^, we have from Eqs. (35), (36), the linearized 

system 

X = y" (49) 

y" = -ax". (50) 

Arguing in the same manner as with Eqs. (37) and (38), wo see that we have 
here a reverse situation, viz., the point Q\ :(— ^,0) is a saddle for a = —k < 
and a center for a = +k > 0. We shall show below that this is indeed the 
case. For this, we shall investigate stability by posing the autonomous system 
as a geometrical Hamiltonian system. The latter technique is said to be more 
reliable than the linearization technique [1] . 

4. Hamiltonian system 

Dropping primes in the autonomous system, Eqs. (35), (36), the Hamiltonian 
system can be defined as 



dH 

dx 
dH 

dy 



= -Y{x,y) = -{ax + cx"^) (51) 
= X{x,y) = y (52) 



The necessary and sufficient condition for the system (51), (52) to be a Hamilto- 
nian system, namely, ^ + ^ = 0, is fulfilled for all x and y. [Such fulfillment 

is a special feature of the GR path Eq.(4)]. Moreover, ^ = and therefore 
H{x,y) =const. independent of <p. Prom the Eqs. (51) and (52), we get 

Hix,y) = ~x''-^x^ + uiy) (53) 

H{x,y) = ly'+<x) (54) 

where u{y) and v{x) are arbitrary functions subject to the consistency of Eqs. (53) 
and (54). These two equations will match only if 

u{y) = \y^-C (55) 
v{x) = -^x'-^x'-E (56) 

where E is an arbitrary constant. The Hamiltonian paths are given by 

H{x,y) = -^x^-^x'+^-y^-E (57) 

where E is a parameter. It follows that 
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dxdy 



= -{a + 2cx) (58) 
= 1 (59) 
= 0. (60) 



As before, the equilibrium points occur when X = and y = which give the 
points Pi : (0,0) and Qi : (—-,0). Thus the quantity 



d'^H d^H 



dx"^ dy"^ \ dxdy J 



has the following values 



50 I Pi = -01 (62) 
?o I Qi = a. (63) 

When —1 < Q: < 0. the equilibrium point Pi is a stable center since qo > 0, 
but Qi is an unstable saddle point. For a > 0, the conclusions are reversed. 
These confirm the results of Sec. 3. The value a = is a bifurcation point as it 
represents a transition of the system from a stable center to an unstable saddle 
and conversely. 

We shall now see what result do we get applying the present method to light 
trajectories for which a = 0. From the original set of Eqs. (15), (16), we see 
that they lead to the same Hamiltonian set of Eqs. (51), (52) with the difference 
that a is now to be replaced by h. The equilibrium points then are P2 :(0,0) 
and Q2 :(— 1,0). The point P2 implies 

a; = r ^ 00 (64) 

and the value of go |p2= —b= 1 > 0. This implies that P2 is a center. From 
this, we learn that light trajectories (straight lines) in asymptotically flat space 
(r — » 00) are stable. This is an expected result. The other equilibrium point Q2 
implies 

x = --=>r = 3M (65) 

c 

at which qo \q2— b — —\ < Q showing that Q2 represents a saddle point. In 
other words, light orbit at r = 3M is unstable. This shows that the instability 
of circular orbits of light at i? = 3M depends only on the sign of h and is 
independent of the sign of a, unlike in the case of material orbits. 

5. Case 11: X^Q 
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Wc have to consider the full autonomous system (11), (12) and as usual, the 
equilibrium points are given by x = 0, y = 0- The latter gives the equation 



g{x) = cx^ + bx'^ + ax'^ + d=Q,xj^ 0. (66) 

Since x has to be non-negative, we have to look only for positive roots of ^(x) = 
0. Once the known signs of coefficients are plugged into g{x) = 0, we may 
apply Descartes' rule of signs to see that g{x) =0 can have either one or three 
positive real roots, the rest are either negative or imaginary. The auxiliary 
equation ^'(a;) = ^ = has two zero roots and two nonzero roots /x^, given 

by 

-26 - V46^ - 15ac -26 + \/46^ - 15ac , , 

Mi = ,M2 = • (67) 

The reality of the roots of g'{x) = demands that 7 = 46-^ — 15ac > 0. Since 
6 = — 1 and c > 0, we see that /x^ < /i2. Let us denote a representative positive 
root of g{x) by 77 7^ 0, that is, 9(77) = 0. The representative equilibrium point 
is then {x,y) = (r?, 0). Then we employ the usual operations on Eqs.(ll) and 
(12), viz., a translation x' = x — r], y' = y, followed by linearization in the 
neighborhood of {x',y') = (0,0). The final result is 

x = y' (68) 

y = y ^ y. (69) 

Using Eqs.(39)-(42), we get 

^ ^ _ (W±j^l±^^ ,p = 0, A = 4 (5£^!±i^!Lt^) . (70) 

For a meaningful analysis, we must have 7^ which means g'ij]) ^ 0, that 
is, 1) can not be a repeated root of g{x) = 0. Thus, we find that {x',y') = 
(0, 0) will be a saddle if g < and A > 0. This is possible if either r] < fi^ or 
Tj > IJ.2- The point (x',y') = (0,0) will be a center if 17 > and A < which 
means /U^ < 77 < /Zj. The linearization scheme is not applicable for 77 = 0. The 
important point to note here is that /Lti,/X2 do not depend on the cosmological 
constant A. Thus the constraint 7^ = 4 — > applies to orbits resulting 

from the effect of M alone. Orbits close to the static radius are not sensitive to 
this constraint. 

With the above general picture in mind, let us numerically study the be- 
havior of approximate roots of g{x) = for some choices of h. For a given 
M the equilibrium points, hence the radii, vary depending on the values of h, 

Xeq{h,X). Choosing units in which M = 1, and with the values of 

A 

3h^ 

can be rewritten as 



coefficients given by c = 3,6 = — l,a = h = — grr the equation ^(a;) = 



h'^{3x^ -x'^)+x^ - ^ =0;h^0 (71) 
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We observe the following behavior. When — > 0, we get only one very small 

root that approximates to the static radius Xeq = Xgt = ) ^ ~ 10~^*. Other 
roots are imaginary. As we increase h? up to ^ , we see that the picture remains 
almost the same, that is, wo continue to obtain a single radius of the order of 
Xst- When /i^ > ^ or 7^ > 0, we find that there occur three positive roots, one 
is of the order of the same static radius, but the other two roots correspond to 
orbits in the vicinity of M. These results confirm that the radii of orbits close 
to or at the static radius are indeed insensitive to values of 7. Let us consider 
a specific value h = 8 (say), then we have the following equilibrium points: 
Pi : {x = Xst, y = 0), P2 ■■ [x ^ 0.016, y = 0) and P3 : (x = 0.316, y = 0) while 
III = 0.012, = 0.254. According to the general discussion above, we expect 
that P2 should be a stable center as < x < ^2 while Pi and P3 should be 
unstable saddles. 

Let us confirm the results by the method of Hamiltonian system. Following 
the same procedure as in Sec. 4, we deduce that 

H{x,y) = ly' - {ax + + \x^ - \x-^) - F. (72) 

where F is an arbitrary parameter. The expression for go is 

go = -{h^lcx-Mx-"^) (73) 
= 1 - 6a: - ^x-^. (74) 

Prom this, we can conclude the following: At the lower bound, that is, at the 

local photon radius Xph = 1/3, wc find qo — —1 since h'^ — 00. Therefore 
this particular orbit is unstable and the instability is independent of A. On the 
other hand, at the static radius, go \x=Xst> 0) implying that the photon orbit 
(again h"^ = 00) is stable at a; = Xgt- The stability of light orbits at the static 
hypersurface is similar to that in the asymptotically flat region discussed in 
Sec. 4. At X = go < 0, hence P = 6 is also an unstable radius. Furthermore, 
go |p2> and go |p3< confirming earlier expectations. At the static radius 
X = Xst, we have h"^ = 0, and Ax^** ~ 10^^ so that go \pi= —00. This shows 
that circular material orbit at the static radius (Pi) is unstable. 

What then is the upper bound Rub for stable circular material orbits? This 
can be found by requiring that go > or 

6a; + ^x-^ < 1. (75) 

Putting the expression for from Eq.(14), and assuming that Rub '> 6, we 
find that 
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Rub = 4-3 j (76) 

which is slightly smaller than rgt- The radii at which orbits begin to be stable 
can be obtained from go = which gives 

A = /i2a;^(l-6a;),a;7^0 (77) 
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The maximum of A is located at a; = ^ . Again using the expression for h?' from 
Eq.(14), we get 

Amax = ^ - 0.000711. (78) 

For A ^ 10^^'"' cm^^, the maximum Schwarzschild mass is Mmax = (^^a-)^ ~ 
5.75 X 10^°M©. Stable material circular orbits can exist only at i? > ^ cor- 
responding to 7 > 0. This is confirmed by the stability at P2 {R — 10) and 
instability at P-i{R ~ 3.3). The values A,„ax and R = ^ correspond to another 
critical value h"^ = as may be obtained from Eq.(14). We have obtained it 
here from a totally different consideration, namely, of roots of ^(x) = 0. For 

> X' there exist local stable equilibrium orbits x^q while no stable local 
Xgq exist for < We see that the restriction /i^ > ^ is weaker than the 
previous >V2 for the case A = 0. 

It is remarkable that a single parameter (70 completely reproduces all the 
results obtained by Stuchlik and Hledfk [3], including their numerical value, viz, 
Vcims) = ^ = 0.000237. 

The phase space method can also be applied in the pure de Sitter space 
which corresponds to A > 0, M = 0. There is now no balance of forces at any 
radius, hence there is no static radius. The metric with M = immediately 
fixes a = c = in g{x) = 0. The equilibrium points then occur at x'^ = — 
This implies that there are no real equilibrium points and we conclude that 
circular orbits are not possible in this space. 

6. Path equation in Kerr spacetime 

The phase space method can be profitably utilized in the study of motions 
in a refractive medium as well. For instance, Evans and Rosenquist [4] showed 
that the equation of optics in a refractive medium of index n{r ) can be effec- 
tively rephrased as a Newtonian " / = rna" form of mechanics. This optical- 
mechanical analogy led via Fermat's principle to a path equation for light in the 
form 

where ^ is a stepping parameter defined in Ref.[4] by dA = n~'^dt. Optical 
analogues of mechanical quantities are marked by The equation of motion 
(79) has been subsequently extended in Ref.[5] to include also the motion of 
material particles. Note that the form of n{l^) can be arbitrarily prcassigned 
depending on the nature of the medium. The relevant quantities in this formal- 
ism are the optical version of mechanical quantities. For instance, instead of 
the classical angular momentum h, its optical analogue "ho'\ viz., ho = r'^^ 
is conserved if n = n{r). A specific form of n = n(r) depicting a Schwarzschild 
gravitational "medium" exactly yielded the path Eq.(4) for A = [6]. (See also 
Ref.[7] for another interesting derivation). It is clear that complicated forms of 
n{r) corresponding to arbitrary spherical media would lead to path equations 
more complicated than Eq.(4). In these cases, the present method might be 
preferable to conventional methods. 
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An example is Kerr spacctime which represents a unique rotating black hole 
solution for A = 0. Alsing [8] has extended the "medium" analogy to Kerr 
spacetime with rotation parameter J (= angular momentum per unit mass M 
of the rotating source) and obtained, to first order in — , the following path 
equations on the equatorial slice: 



■u- 3Mu^ = — r 



J.2 



1_ - I ( i_8Mu^V2^ 



{part(^0) 



= -^{light) (81) 

Here vq is the initial velocity of the particle at infinity and L is its conserved 
total "angular momentum" per unit test mass given by 

^-^dA-^ ^^^^ 

where p = re"^*". (In the asymptotic region, n = 1, A = t p = r so that, 
for J = 0, we have L = h, the familiar mechanical angular momentum). For 
a particle starting a radial fall from infinity, L = 0. In this case, with n{r) ~ 
1 + 2M g obtains to lowest order in ^ that ^ = This imphes that 

p ' r at ^ 

the particle starting with an initial radial fall begins to co-rotate with the black 
hole in its vicinity (Lense-Thirring effect). In general, we shall take L ^ 0. 
Conventionally, the sign of J is taken as positive or negative according as the 
source rotation is in the counterclockwise or clockwise sense. 

To examine stability, we first note that the autonomous system is of the 
same type as in Eqs.(15), (16), only the coefficients are different. The next step 
is to follow the same procedure as in Sec. 3. Applying it for a particle starting 
at vo = 0, we get the equilibrium points in the {x',y') plane at Qi : (0,0) and 

Q2 : (-f,0) where 



/? = ±\/(l + -^j --^ (83) 

which reduces to a when J = 0. The interesting result is that the reality of 
P immediately imposes two restrictions, viz., that J ^ ^^'^ 
quantity under the radical sign in Eq.(83) must be positive which implies ( J — 
J+){J — J-) > =J> either J < J± m J > J±. These restrictions must be 
respected if circular orbits are to exist at all in Kerr spacetime. Once this is 
fulfilled, exactly the same arguments about the stability as in Sec. 3 go through 
under the replacement of a by /3. 
Case (a): = 

This degenerate condition corresponds to two critical values J+, J_ of J 
which are 



J±=T 



8M2 



2V3^±1 

1j 



(84) 
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The phase paths arc the same as those given by Eq.(25) indicating unstable 
equiUbrium at radii given by 



Kerr _ /oc^. 



Putting the values of J± from Eq.(84), we find the radii r^^^^ = =F\/3-^, which 
implies that the rotation J of the source has no role in determining the radii of 
circular orbits if 0^ = 0. 
Case (b): 0^ > 0. 

The phase paths are the same as Eq.(36) with the replacement of a by p. 
Hamiltonian analysis reveals that Qi is a center and Q2 is a saddle if /3 < 0. 
These conclusions are reversed if /3 > 0. Stable circular orbits occur at 



r 



err 



6M 



L2 



(86) 



Prom the above, it follows that, as J ^ ±00, the stable radii r^_|g'' go far beyond 
6M. This implies that the rotation J of the source can not bring about stable 
circular orbits at radii below 6M for material test particles. The situation is 
the same as in the nonrotating case. 

Case (c): a ^ 0. 

For light (fo = co), a = ^' hence the equilibrium points occur at 

Ri : (0,0) and i?2 : (-^,0) where 



, /. 24JM2 

Arguments similar to Case (a), Sec. 3 go through, with the replacement of a by 
(7. Thus stable radii occur at 



6M 



r^rjkt = ^— (88) 



depending on the sign of a. When J = 0, the orbit has infinite radius. When 
-1 < < 0, stable circular orbits will exist for r^fJ^f[^ > 6M but when 

^"^^ > 0, there can not be any stable radius for light because r^ig^t becomes 
negative. 

7. Optical black holes 

The advantage of the phase space method is that it can be applied to situa- 
tions beyond known gravitation theory when a potential function is not always 
evident. This can occur, for instance, when one deals with a path equation in 
the environment of a simulated black hole described by a refractive medium 
with index n(r). Possibility of laboratory creation of such optical black holes 
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exist in view of a remarkable experiment [9] performed in Bose-Einstein con- 
densates. The experiment demonstrated that optical pulses can travel in the 
condensate with extremely small group velocities, as low as 17m/ s. The group 
velocity of light in the vicinity of a real gravitational black hole can indeed be 
arbitrarily low [10]. It has been shown that light motion around a dielectric 
vortex structure mimics motion around a black hole [11]. Creation of an event 
horizon would require that the vortex flow be supplemented with a radial flow 
as well [12]. Interesting physical effects, like optical Aharonov-Bohm effect far 
away from the vortex core and bending of light near the core, stem from the 
consideration of a dielectric medium having a velocity field 'u and a varying 
index of refraction n. 

One might obtain the trajectory of light directly from Fermat's principle 

8 j n{'r)dl = 0. (89) 

The resulting path equation is given by Eq.(79) which, on the equatorial plane 
9 = 7r/2, gives 

dr [r^'n?{r) — r'lr'^]'^ 

n 

where is a constant of integration, interpreted as impact parameter. Stability 
of circular orbits can be easily studied directly once a form for n(r) is given. 

A plausible form for n(r) simulating a static dielectric medium of optical 
black holes has been studied by Marklund, Anderson, Cattani, Lisak and Lund- 
gren [13] . It is given by 

n2(r) = l + ^ (91) 

which has a divergence at r = 0, and is a constant. For r 3> rg, n(r) ^ 1 and 
for r <C ro, n(r) ~ ^. Defining, as before, u = x = ^, the autonomous system 
corresponding to the problem can be written as 

X = y (92) 



y = ex,e= ^° ■ (93) 
rf 



The equilibrium point (0, 0) refers to circular orbit only at the asymptotic re- 
gion r = 00. Taking the limit oo, we find e = —1. Hence the path 
equation — ex^ = C represents a family of concentric ellipses around the 
origin showing that the orbit is stable independent of any finite value of the 
extent ro of inhomogeneity. But no circular orbit at a finite radius is possible. 
This result is quite consistent with the nature of various trajectories analyzed in 
Ref.[13]. However, a pathological solution of y = may be imagined by taking 
e = Tq = r^. We find that there is no equilibrium point at all in this case 
as i ^ although y = 0. The latter yields a phase path equation y = ^ = \fC 
which integrates to give real space trajectory x cc (j). This is just in the form of 
Archimedes' spiral Vi/r = (j) as shown in Ref. [13] . 
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Wc show now that circular orbits arc possible in a nonuniformly moving 
medium with a slowly varying refractive index. Under these conditions, Leon- 
hardt and Piwnicki [11] considered a vortex core with a velocity profile decaying 
away from the core 

^ W 

" = (94) 

where 2'kW is the vorticity. Let us formally introduce the index given by Eq.(91) 
into the Hamilton- Jacobi equation for light motion derived in [11]. The resulting 
path equation can be translated to the following autonomous system in the far 
field limit 

i = y (95) 
y = -a: + -3!- (rgx + 2r^H^V + 3r^WV) (96) 

''AB 

where Iab is the Aharonov-Bohm modified angular momentum given by 

Iab = 1 + (n^ -l)W = l + -^W. (97) 



r 



Interestingly, we find that there is only one equilibrium point (apart from the 
trivial one at a; = = or r = 00) at a finite radius photon orbit 



R=\/3ro 



(98) 



In order that this radius be real, we must have 

/2 _ 



1 I 'WAB-rl) 



W2 -^'>1 (99) 



where A'' is real. The quantity qo in this case works out to 

4\ 



qo = - (^3 j - > (100) 

if > 1. Thus the orbit is stable. Of course, the conclusion crucially depends 
on the form of n{r). 

We finally mention an interesting similarity between L of Eq.(82) and Iab 
of Eq.(97). A light ray approaching radially (l — 0) will acquire an angular 
momentum (n^ — 1)W near the vortex core very similar to the Lense-Thirring 
effect. With ~ 1 + where m is some constant, we have ^ ~ ^^J^^- This 
allows us to identify 2W as the angular momentum of the vortex motion. 

8. Summary 

Phase space analysis has been successfully applied to practically all walks 
of life, from physics, engineering, biology to social sciences. Somehow its use in 
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gravitational physics seems rather scarce. Our motivation here was to fill that 
gap. There is certainly room for further development like exploring how other 
sophisticated techniques from the phase space repertoire could be applied to 
stability of noncircular orbits or even classical fields. 

We obtained information on stability of circular orbits arguing from the 
geometrical path equation alone. The usefulness of the method is demonstrated 
in several situations of physical interest. When the cosmological constant A = 0, 
the dimensionless parameter a played a key role in determining the stability of 
actual orbits in the physical (r, (f)) space. It was shown that a = is a cusp 
describing a marginal state, viz., the radius R = 6M is neither stable nor 
unstable. The light orbit at i? = 3M is unstable independent of the sign of 
a and in the asymptotic region it is always stable. All the conclusions were 
confirmed by the method of Hamiltonian system. 

When A ^ 0, there occurred either one or three equilibrium points. The one 
corresponding to static radius does not depend on the parameter 7 but depends 
on A. We found that, for A > 0, trajectories of circular material orbit at the 
static radius are unstable. However, for light orbits, h 00, so that go = 1) 
hence circular light orbits at static radius are always stable. These two results 
help us understand better the nature of the static hypersurface: Even though 
forces balance at the hypersurface, it is not exactly like the usual flat asymptotic 
region where both matter and light orbits are stable. The other two equilib- 
rium points arc local and already analyzed in Sec. 4. They do not depend on A 
implying that circular orbits in the vicinity of M are not influenced by A. This 
is a physically consistent result. It was shown how a single parameter qo nicely 
reproduced all the relevant results about circular orbits in the Schwarzschild-de 
Sitter spacetime. 

We dealt with circular motion of light and massive particles in the equatorial 

plane of the Kerr black hole. The derived results are new. The restrictions 
J < J± OY J > J± for the existence of circular orbits were derived. A general 
parameter /3 was found that reduces to a under zero rotation, J = 0. For 
the value = 0, the rotation J of the source has no role in determining the 
radii of circular orbits. For > 0, we found that even high values of source 
rotation can not bring about a stable radius below 6M. These results could be 
of importance to accretion phenomenon in astrophysics. 

Finally, we applied the method to real optical dielectric (static and mov- 
ing) and obtained consistent results. Here we only discussed a simple example 
relevant to optical black holes but any given form of refractive index n(r) can 
be similarly handled. The important advantage is that we did not require in- 
formation on dynamical potential functions, but relied solely on the path equa- 
tions coming from Fermat's principle or Hamilton- Jacobi equation. Several such 
path equations corresponding to various refractive media have been worked out 
in Ref . [4] . It will naturally be of interest to apply the method in the refractive 
wormhole "media" constructed from exotic matter [14] or in Brans-Dicke theory 
[15,16]. Work is underway. 

Figure captions 

Fig.l. The origin O is a cusp [See p. 19, Ref.l]. Path AO leads to the origin 
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while OB leads away from O. These paths correspond to C = 0. For C 7^ 0, the 
paths CDE never reach O and the motion is analogous to the whirling motion 
of the bob of a pendulum. 

Fig. 2. The origin O is a stable center. Initial conditions slightly shifted 
from the center take the phase paths on closed elliptic orbits around O. The 
corresponding motion in physical space is periodic analogous to small oscillations 
of the bob about downward vertical. 

Fig. 3. The origin O is a saddle point. Only the paths AO and BO approach 
the origin while OC and OD move away from it. Other paths do not lead to 
the origin. EF, GH represent whirling motion. 
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